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ABSTRACT 

It has long been suspected that metal polluted white dwarfs (types DAZ, 
DBZ, and DZ) and white dwarfs with dusty disks possess planetary systems, but 
a specific physical mechanism by which planetesimals are perturbed close to a 
white dwarf has not yet been fully posited. In this paper we demonstrate that 
mass loss from a central star during post main sequence evolution can sweep 
planetesimals into interior mean motion resonances with a single giant planet. 
These planetesimals are slowly removed through chaotic excursions of eccentric- 
ity that in time create radial orbits capable of tidally disrupting the planetesimal. 
Numerical iV-body simulations of the Solar System show that a sufficient number 
of planetesimals are perturbed to explain white dwarfs with both dust and metal 
pollution, provided other white dwarfs have more massive relic asteroid belts. 
Our scenario requires only one Jupiter-sized planet and a sufficient number of 
asteroids near its 2:1 interior mean motion resonance. Finally, we show that once 
a planetesimal is perturbed into a tidal crossing orbit, it will become disrupted 
after the first pass of the white dwarf, where a highly eccentric stream of debris 
forms the main reservoir for dust producing collisions. These simulations, in con- 
cert with observations of white dwarfs, place interesting limits on the frequency 
of planetary systems around main sequence stars, the frequency of planetesimal 
belts, and the probability that dust may obscure future terrestrial planet finding 
missions. 

Subject headings: planetary systems-white dwarfs-numerical techniques: N-body- 
asteroids 
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Introduction 



Hundreds of exoplanets have been discovered around a variety of stars, from small M 
dwarfs out to large G giant stars, to tight orbits around pulsars. The conclusion is clear that 
planet formation is a robust process in our galaxy. A corollary to this conclusion is the idea 
that observable evidence for planetary systems must therefore exist from early star formation 
to various stellar end states. The challenge is to identify which phenomena are uniquely due 
to planets, and which are considered false positives. When planets are discovered at each 
stage of stellar evolution, it will be additionally challenging to interpret the impact stellar 
evolution had on the surfaces, composition, and orbits of such planets. Studying these 
physical processes tell us what the eventual fate of our own Solar System will be. 

A very compelling case can be made for the presence of planetary systems around metal 
enriched white dwarfs and white dwarfs that possess strong infrared excesses. Many studies 
with Spitzer and ground based optical spectroscopy have provided nearly twenty dusty white 
dwarfs with some sort of abundance estimate for the dusty material, either through detection 
of the 10/im silicate feature, careful abundance analysis of narrow metal lines in the pho 
tosphere of the target w" 
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2011bJ). More recently, the Wide Field Infrared Survey 
201 Id ), has shown that ~1% of WDs show dusty disks, con- 



firming earlier work with Spitzer. In contra st, nearly 20-25% of all WDs show measurable 



metal accretion (e.g. IZuckerman et al.ll2003l ). 



There have been two main conclusions from this excellent body of work: 1) that the 
material drizzling onto the surfaces of these white dwarfs is primarily similar to the inner 
terrestrial planets in our solar system, and 2) that planetesimals wh ich are tidally disrupted 
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Jural (120031 ) and lJura et al.l (120071 ) have elegantly shown that most white dwarfs with 
infrared excesses require inner disk radii of just a dozen or so white dwarf radii, or about 
0.1 R . The extreme proximity to the white dwarf produces quite unusual conditions for a 
dust disk, and argues strongly for a tidally disrupted asteroid that has arrived from beyond 
a few AU, since anything at smaller distances during post main sequence evolution would 
be evaporated within the envelope of the evolving giant star. Based on these arguments, it 
has generally been assumed that a planet several AU away from the WD is responsible for 
perturbing planetesimals into highly elliptical, tidally disrupting orbits. 



A primary uncertainty in linking planetary systems to dusty and polluted white dwarfs 
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is the exact d ynam ical mechanism that delivers planetesimals so close to the white dwarf. 



Alcock et al.l (119861 ). in a prescient paper, tried to explain the presence of the metal polluted 



hydrogen WD (type DAZ) G74-7 with infrequent cometary impacts onto a white dwarf 
surface from Oort cloud analogues. Further work was done looking at the dynamical stability 
of planetary systems during post-main sequence evolution as a possible driver for material 
into the inner system (IDebes fc Sigurdssonll2002l ). This scenario could reproduce the rough 
number of polluted white dwarfs assuming that as many as 50% of planetary systems are 
unstable. It also predicted that some planetary systems may be stable for up to 1 Gyr 
after a star turns into a WD, in accordance with the range of pollution observed for WDs. 
Because dusty white dwarfs are observed at relatively late cooling ages when the WD has 
very low luminosity, processes that are initially efficient at creating dust durin g post main 



sequence evolution become n egligible on ce a white dwarf is fo rmed and is cooling (IDong et al. 



2010l ; iBonsor & Wyattll20!0l ). Recently, iBonsor et al.l ( 1201 ll ) investigated the frequency with 



which planetesimals in belts exterior to a giant planet might be perturbed into the inner 
system. They found in their iV-body simulations that a sufficient number of planetesimals 
can be perturbed to match the observed evolution of the dust accretion rate onto white 
dwarfs as a function of cooling age. Two main assumptions however were used in that 
work: that the planetesimals must survive post main sequence evolution, and that roughly 
10% of perturbed planetesimals would be further deflected (presumably by a second interior 
planet) into tidally disrupting orbits. Most planetesimals beyond the ice line (Kuiper belt 
type analogues) should be primarily icy. They may significantly sublimate during post main 
sequence evolution before they can be perturbed. Similarly, many interior planets that could 
further perturb outer planetesimals into white dwarf-crossing orbits may be engulfed during 
post main sequence evolution. Both scenarios could be helped by eccentricity pumping of 
icy bodies further away from the i nner system that don't see the mass loss of the central 
star as adiabatic ( iVeras et al.l 120 111 ). 

While all dusty white dwarfs show evidence for metal accretion onto their surfaces, 
not all metal enriched white dwarfs show evidence of a dusty disk. A correlation between 
accretion rate and the presence of dusty disks has been suggested (Ivon Hippel et al.l 120071 ; 
Farihi et al.l l2009t IZuckerman et al.ll2010f ). with an accretion rate of >10 8 5 g/s consistent 



with the range above which dusty white dw arfs occur. The delineation could be due to the 
average size of the perturbed planetesimal ( jJural 120081 ). In this scenario, dusty disks only 
occur for larger than average disruptions, while smaller planetesimals that are disrupted are 
quickly turned to gas through mutual collisions and sputtering. Evidence for the presence 
of multiple smaller disruptions causing gaseous disks h as been shown with the discovery of 
circumstellar Ca gas absorption around WD 1124-293 ( IDebes et al.l 1201 lai ). 



In this paper we present a new scenario for perturbing planetesimals into highly eccen- 
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tric orbits, thus strengthening the link between dusty white dwarfs and planetary systems. 
We suggest in §2] that perturbations in eccentricity from interior mean motion resonances 
(IMMRs) with a giant planet roughly the mass of Jupiter are sufficient to create a steady 
stream of white dwarf crossing planetesimals. In particular, we hypothesize that the 2:1 res- 
onance is most efficient at driving white dwarf crossers, and that these asteroids are quickly 
tidally disrupted. We use numerical simulations to follow the dynamics of the Solar Sys- 
tem's asteroid belt under the influence of post-mai n sequence evolution and perturbation by 
Jupiter in §21 Our model differs from the results of iBonsor et al.l (120111 ) primarily in that we 
choose p lanetesimal s that should be primarily rocky and survive post main sequence evolu- 
tion (see I Jural 120081 ) . and we follow the dynamics of the planetesimals from perturbation to 



entering the tidal disruption radius of the white dwarf. We also simulate the close approach 
of an asteroid to a white dwarf in §3] to determine whether asteroids can be tidally disrupted 
quickly. In §H] we confirm that a significant fraction of asteroids are perturbed into close 
encounters with a white dwarf, and that highly eccentric encounters between a small rubble 
pile asteroid and a white dwarf are sufficient to tidally disrupt the asteroid. We compare 
the results of our models to the currently known population of dusty and polluted WDs in 
§3 including estimating the asteroid belt masses necessary to produce the observed metal 
pollution in white dwarfs, and discuss our conclusions in §61 



2. The Interior Mean Motion Perturbation Model 



The Kirkwood gaps of the Solar System's asteroid belt are regions where asteroids 
are quickly removed due to gravitational interactions with planets. Within IMMRs, an 
asteroid's eccentricity random walks until it undergoes a close encounter with a planet (or 
planets) and is ejected, col lides with a planet, or collides with the central star ( jMorbidelli 



19961 ; iGladman et al.lll997l ). This motion is limited to those bodies in an IMMR, whose 



width 5a max , can be roughly approximated by the maximum libration width for interior first 
order resonances in t he restricted three body c ase and expanding the equations of motion 
for low eccentricities (jMurray fc Dermott 2000 ): 
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5a max = ± ( y— e ) ( 1 + 07 . 2 3 - 7T. -a, (1) 

where n is the mean motion of the planetesimal and e is the eccentricity of the planetesimal, 
C r is a constant from the resonant part of the disturbing function. The constant comes 
from the resonant argument and determines which resonance is being used for a calculation, 
and a is the ratio of the asteroid's semi-major axis to Jupiter's. The quantity — is given 
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by: 



C 

— = t*a\fd(<x) 
n 



(2) 



For the 2:1 resonance, and afd (a) =-0.749964 (IMurray fc Dermottll2000i ). As can 

be seen in Equation (H the maximum libration width implicitly depends on the mass ratio 
fi. As a star evolves and loses mass, the mass ratio /i increases. As a result, the width of the 
resonance 5a max increases and bodies previously exterior to t he IMMR become trapped. Just 



as in Hill stability and the stabil ity of multi-planet systems ( )Gladmanlll993l ; IChambers et al. 



19961 ; iDebes fc Sigurdssonll2002l ). mass loss from the central star increases the perturbative 
influence of the planet. Figure [1] demonstrates the growth of Sa max with a corresponding 
change in [i for the 2:1 resonance. 

The 2:1 resonance in the Solar System appears to be the most efficient mechanism for 
scattering asteroids into Sun crossing orbits over Myr timescales. iV-body simulations of 
asteroids injected into different resonances with Jupiter show a typical lifetime of a few 
million yea rs (such as with the 3 :1 resonance) to a few tens of millions of years for the 2:1 
resonance (IGladman et al.lll997l ). The 2:1 resonance is particularly attractive as a source 
for tidally disrupting asteroids around WDs because the timescale for asteroid removal is 
long (50% removal at >100 Myr), and a significant percentage (~6.5%) of objects in the 2:1 
resonance are perturbed into Sun-crossing orbits. 



We now build on the tidal disruption scenario as envi sioned by lJural (120031 ) . which we 
call the Interior Mean Motion Perturbation (IMMP) model ( IDebesll201ll ). Planetary systems 
with one or more dominant giant planets possess planetesimal belts that are the remains of 
core accretion formation of terrestrial planets/giant planet embryos. Over time these belts 
dynamically evolve through mean motion and secular resonances and lose mass to collisions, 
ejections, and star grazing orbits. Many of the resonances are cleared out to ~ Sa max . During 
post main sequence evolution, smaller planetesimals migrate and evaporate in response to 
gas drag from the central star evolving off the main sequenc e, but a sig nificant portion of 
large (r>l-3 km) planetesimals survive primarily unscathed ( 1 Jural 120081 ) . Meanwhile, due 
to mass loss from the central star, the reservoir of planetesimals that are now vulnerable to 
eccentricity excursions from mean motion resonance pert urbations has grown, either through 



trapp ing of medium sized planetesimals in resonance ( IDong et al.l I201CH ; iBonsor &: Wyatt 
20101 ) or through the growth of <5a max - 



Planetesimals that are on initially circular orbits interact gravitationally with a domi- 
nant giant planet and become perturbed into a highly elliptical orbit. The perturbed plan- 
etesimals eventually obtain such high eccentricity that they are tidally disrupted by the 
central WD-the streams of debris from the initial disruption are primarily on similar orbits 
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as their host planetesimal. Over subsequent passes, more disruptions occur, spreading the 
material so that collisions increase. Over many orbits these collisions damp down mutual 
eccentricity and form a coherent disk structure that evolves through a combination of further 
collisions and Poynting- Robertson drag until the dust becomes optically thick. At this point 
the disk evolves viscously until it becomes optically thin, or a ll dust is s puttered into gas 
from grain-grain collisions and the debris from other asteroids (I Jural 120081 ) . 



3. Numerical Methods 

In order to test the above IMMP model, we will demonstrate that planetesimals 1) 
entered a highly eccentric orbit within the tidal disruption radius, and 2) disrupted within 
the expected disruption radius. Ideally our simulations would also follow the disruption at 
later times to see if a dust disk forms that matches observations of known dusty white dwarfs, 
but this is computationally intensive with the techniques we use in this study and is saved 
for future work. 

To demonstrate that the IMMP model satisfies the first criterion, we dynamically sim- 
ulated the response of the Solar System's asteroid belt to post main sequence evolution. 
We chose the Solar System primarily because it should be a decent proxy for a planetary 
system where only one planet dominates the inner planetesimal region, and should not be 
subject to any biases associated with an incorrect treatment of Gyr of dynamical evolution 
between a planet and its planetesimal belt. Our simulations follow the largest asteroids with 
well known radii-the lower limit to our population is larger than the minimum post main 
sequence survival radius (~10 km) for distances of a few AU. We do not follow smaller aster- 
oids that may be trapped in IMMR resonances during the late stages of post main sequence 
evolution. If we assume that WDs have similar mass asteroid belts, our simulations should 
represent a lower limit to the amount of material tidally disrupted and accreted onto WDs. 

We also simulated the tidal disruption of a small asteroid by a white dwarf using the 
pkdgrav code. We modified the code to treat the approach of a small asteroid to within 
<1R© at various close approaches to determine whether criterion 2) was satisfied, drawing 
on the orbital elements inferred from our iV-body simulations. 



3.1. The Solar System as a Laboratory for Dusty White Dwarfs 



In order to follow the behavior of the 2:1 resonance, we performed numerical simulations 
of large Solar System asteroids in which mass loss from the central star is included. We 
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have take n ten MERCURY simulations using a Burlirsch-Stoer integrator with an adaptive 
timestep ([Chambers! Il999l ) of 710 Solar System asteroids including Jupiter to determine if 
any follow the IMMP model. 

The asteroids were chosen to have radii R>50 km at all orbital semi-major axes and 
R> 10 km with perihelia >3 AU based on the latest known asteroid data compiled by E. 
Bowelfl and the assumption that the Su n will sublimate asteroids smaller and closer than 
this during post main sequence evolution OSchroder fc Connon Smithll2008l ; I Jural 1200 81 ). The 
mass of the Sun was slowly removed over 1000 yr using the equation 



M (t) = 1.0-0.46 



'stop 



--stop 



(3) 



reaching a mass of 0.54 M ( jSchroder fc Connon Smith! 120081 ) . Asteroids were removed 
if they strayed within 1 solar radius and considered tidally disrupted. 

Ten simulations were run for 100 Myr. Four further simulations were performed for 
200 Myr and three separate simulations u sing a Bulirsch-Stoer int egrator from the hybrid 
symplectic integration package described in lStark fc Kuchnerl (120081 ) with a nominal step size 
of l/40th of an orbit of Jupiter and adaptive timesteps were run for 1 Gyr to test whether the 
perturbation declines over timescales comparable to the cooling time of most white dwarfs. 
The 1 Gyr simulations were also used as an independent test to confirm the general behavior 
seen with the MERCURY code. The timescales we cover correspond to a large fraction of 
the cooling age for observed dusty and polluted white dwarfs. For example, a log g=8 white 
dwarf has an effective temperature of ~ 18500 K at a cooling age of 100 Myr, ~ 15000 K at 
200 Myr, and ~8200 K at 1 Gyr. For comparison, all but one of the known dusty WDs have 
T eff >8200 K and 67% of Spitzer observed DAZs have T e// >8200 K. 



3.2. Modeling the Tidal Disruption of a Small Planetesimal 



The tidal disruption is simulated using a "rubble pile" model, based on a collection of 
5000 hard, spherical particles. The simulations use the iV-body code pkdgrav, originally used 
for cosmological simulations, to provide a very fast, parallelized tree for computing inter- 
particle gravity ( iRichardson et al.l 120001 ) . The rubble pile simulations also take advantage 
of its collision resolution, in which particle collisions are detected and resolved according 
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to a number of adjustable physical parameters, namely the coefficients of restitution (both 
tangential and normal). The body is then held together by its own self-gravity with particle 
collisions preventing a collapse. Previously this code has 
pile" asteroid models, in cluding collisions 



20021 ). tidal disruptions (IRichardson et al. 



planet formation simulations (ILeinhardt fc Richardson 
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2005 



Leinhardt et all 120091 ). Tidal 



disruption simulations have proven insensitive to numerous physical parameters, such as the 
tangential and normal coefficients of restitution, as the rotation imparted by the encounter 



and the " depth" of the encounte r with the planet dominate the outcome (IRichardson et al. 
19981 : IWalsh fc Richardso"n1l2006h . 



The modeled tidal disruption consists of a 5000 particle body, having a single close 
encounter with a point mass of 0.5 M Q . Such a point mass has a Roche radius of 89 Rwd ; 
or about 0.9 R©, 695,000 km. We chose orbits with semi-major axes of 4.77 AU, based on 
the orbital elements from one of the tidally disrupted asteroids in our N-body simulations. 
The semi-major axis is larger than one would expect for a Solar System asteroid primarily 
due to the fact that asteroid orbits will expand in response to the mass loss from the central 
star. The asteroid was modeled with different eccentricities to produce encounters inside the 
Roche limit at 80, 75, 70, 65, and 60 Rwd- The progenitor was not spinning on the initiation 
of the encounter, and was only simulated for a single encounter. The simulation was well 
resolved, with a duration of 25,000 timesteps with each timestep representing le-5 yr/2-yr, or 
about 50 sec. At the end of each simulation the state of the aggregate was analyzed, and 
individual orbital elements for each "clump" were produced, along with the statistics for the 
number and size distribution of fragments. 



4. Results 

4.1. The Solar System 

Out of all the IMMP simulations, roughly 2% of the modeled asteroids strayed within 
the tidal disruption radius, with the majority disrupting within the first few hundred Myr. 
The resulting asteroids that were tidally disrupted are plotted as a function of their initial 
semi-major axis a and eccentricity e in Figure |2j Overplotted are calculations of the libration 
width before and after post-main sequence evolution. 

In the 100 Myr simulations, one asteroid impacted approximately every 10 Myr. From 
100 to 200 Myr, the average time between impacts increased to 71 Myr, implying a near 
power law drop in the impact frequency on the white dwarf as time went on. In the 1 Gyr 
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simulations, only one asteroid per simulation was disrupted beyond 200 Myr. In order to 
compare across simulations, the number of asteroids disrupted as a function of time were 
normalized by the total number of asteroids simulated to the per asteroid frequency and 
binned logarithmically in time. Figure El shows the normalized number of asteroids disrupted 
as a function of time for all the simulations. Taken as a whole, this curve would represent 
the average probability at any given time of an asteroid becoming tidally disrupted for the 
size distribution we use. There is a peak at ~30 Myr, which is consistent with a white dwarf 
T e ff ~24000 K. The hottest dusty white dwarfs are slightly cooler than this value and also 
represent some of the largest accretors of material, implying dust rich disks close to where 
the peak of perturbations occur in our simulations. Beyond this peak, the frequency drops 
roughly as t -0 ' 6 , but this is uncertain due to the small number statistics of our simulations. 



We can compare these results to IMinton &: Malhotral (120101 ). who calculated the dif- 
fusion of solar system asteroids (with D>30 km) out of the asteroid belt due to perturba- 
tions from the giant plan e ts. They used second-o rder mixed variable symplectic mapping 
( jWisdom fc Holmanlll99ll ; ISaha fc Tremaindll992l ) to integrate test particles between Mars 
and Jupiter to within 1 AU an d then MERCURY with its hybrid integrator to perform their 
integrations within 1 AU. For lMinton fc Malhotral ( 120101 ). 18% of those asteroids perturbed 
into the inner Solar System impacted the Sun. For our IMMP simulations, 17% of the aster- 
oids that were perturb ed out of the asteroid belt also strayed to within 1 R , our adopted 
tidal disruption radius. IMinton &: Malhotral (120101 ) also found that between 1-4 Gyr, roughly 
0.4% of asteroids in their simulations impacted the Sun, suggesting that even at later times, 
a significant flux of asteroids can be tidally disrupted. Similarly, the late time dynamical 
evolution of the asteroid belt could be approximated by a power-law decline in the total 
number of asteroids close to t -1 , similar to the decline we see in our own impact frequ ency, 
suggesting that the rate of impacts decline with the total mass of the asteroid belt (jJura 
2008h . 



As noted previously in iBonsor et al.l ( 1201 ll ). one can determine an accretion rate M 
onto WDs based on the fraction of asteroidal belt mass scattered (/si), disrupted (/td), and 
eventually accreted (face), as well as the mean time of disruption < t TD >: 



M, 



metal 



/acc/TD/si^belt 
<t TD > 



(4) 



In reality, the accretion rate onto the white dwarf surface is determined by the accretion 
of the dusty or gaseous disks that are generated from the asteroid tidal disruption. Deter- 
mining the precise evolution of a WD dust/gas disk is beyond the scope of this paper, so 
we assume an average accretion rate from the disk is given by < t TD >. We can determine 
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accretion rate vs. time for our simulations if we can calculate a mass disrupted. We took the 
asteroids tidally disrupted in our simulations normalized by the total number of asteroids 
sampled and calculated their mass using their measured radii, assuming a mean density of 
4 g cm~ 3 . This gives us the total mass of ou r asteroids normalized by the mass of the Solar 
System's asteroid belt (Mbeit = 3.6 x 10 24 g) ( jKrasinsky et al.ll2002l ). At each disruption time 
we determined the average mass mass accretion onto the white dwarf by setting < txD > 
equal to the average time between the previous and next disruption, or the end time of the 
simulati on and assum i ng th at f acc was of order unity. This mostly follows the same proce- 
dure as iBonsor et al.l (120111 ). with the exception that we implicitly determine the quantity 
/TD/krMaeit by following individual disruptions, and we assume that / acc is more efficient. 

In principle the results of our simulations can be scaled to any starting asteroid belt 
mass by the total mass of asteroids, assuming a similar dynamical architecture of a dominant 
giant planet with an interior belt of asteroids following a similar size distribution to that of 
the Solar System's asteroid belt. Figure 0] shows our results for the accretion evolution of 
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Farihi et al. 
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2009. 


2010) 



simulations into corresponding T e ff assuming a hydrogen WD with log g=8.0. Looking at 
all disruptions in our simulations, we take a power-law fit to our results and find that the 
mass accretion rate drops as t~ 2 . We can extrapolate this relationship between time and 
accretion rate to the oldest WDs. This is shown as the solid line in Figure HI The inferred 
accretion rates depend on our assumptions in the following way: 



M metal - M metal)0 (i^j (3 5^4 J • ( 5 ) 

Our results suggest that the required mass in asteroids around other white dwarfs might 
have been much larger than our own Solar System, by as much as a factor of 10 3 . However, 
several uncertainties exist in our model. The biggest uncertainty for this model is the exact 
location of the giant planet and how the timescale for asteroid perturbation out of the 2:1 
resonance might scale with planetary architecture. We also haven't accounted for how the 
accretion rate might evolve over the lifetime of a gaseous or dusty disk caused by the tidal 
disruption of an asteroid. One would expect orders of magnitude higher accretion rates soon 
after a tidal disruption, with a power-law or exponential decay in accretion rate as mass 
is lost from the disk of material that forms from a tidal disruption-this effect would be 
balanced by the fraction of a WD's cooling age over which this high accretion phase might 
occur. Such modeling would require a more careful analysis of how such disks evolve and 
the main mechanism for accretion onto the WD. Recently, it has been shown that Poynting- 
Robertson drag may drive accretion in some cases, but cannot account for all of the accretion 
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observed flRafikovll2011al Jbl: ljura fc Xulboilf ) 



It is also unclear how tidally disrupted streams of material will evolve into dusty /gaseous 
disks and over what timescales this process occurs. If the relative velocities of collisions are 
too high close to the WD, dust grains will suffer evaporative collisions and might not be effi- 
cient at leaving material close enough to accrete onto the WD. Similarly, if relative velocities 
are too high, the tidally disrupted streams of material may not be able to collisionally damp 



into a disk with low enough eccentricity to efficiently accrete ([Shannon fc Wull201ll ). Our 
sample of asteroids encompasses most of the largest asteroids present in our Solar System, 
but is missing most asteroids with radii <20 km. Therefore it is possible that our implied 
tidal disruption rates could be higher if we take these objects into account, or other planetes- 
imal belts have different size distributions such as what might be expected if a l arge number 



of pla netesimals are trapped in resonances during post-main sequence evolution (IDong et al. 



2010). These uncertainties affect many of our assumed model values for the efficiency with 
which tidally disrupted material is accreted onto the WD and the average timescale between 
tidal disruptions. Further work on these uncertainties through better observations of dusty 
disks, the modeling of WD disks, and longer simulations of dynamical perturbations with 
more particles could help to better predict how efficient the IMMP model is at delivering 
material close to a WD. 



4.2. Tidal Disruption 

In all of our tidal disruption simulations, the asteroid was significantly affected by the 
tidal forces of the central white dwarf, with varying degrees of disruption after the first 
pass. As the periastron decreased, more violent disruptions occurred, with more mass being 
generated in smaller bodies. 

Figure [5] shows a snapshot of one asteroid disruption for a close encounter to 60 -Rwd- 
The asteroid is elongated at the early stages of the disruption and quickly spreads out 
over 2 AU from head to tail of the train of debris after the first pass. Despite this, the 
disrupted stream is rather similar in its orbital elements. This tight stream is reflected for 
all the encounters and the orbital elements of the fragments were tightly centered around the 
original orbit of the incoming asteroidal body. Further work needs to be done to determine 
whether these fragments will settle through mutual collisions and further tidal disruptions 
into a more circularized disk. This suggests a phase of evolution, perhaps lasting many 
orbital timescales for the fragments, where dust is in an elliptical distribution. In fact, 
observatio ns of dusty white dwa rfs with emission lines from gaseous disks show evidence of 



ellipticity (IGansicke et al.l 120081 ) . though not at the level implied by our simulations, where 
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the eccentricity of the debris stream is in excess of 0.99. Over many orbital timescales, 
the fragments will spread and precess at slightly different rates and start colliding. The 
exact evolution of the system will depend on the timescale for fragments to be directly 
perturbed onto the WD surface by the giant planet and the timescale for mutual collisions 
and a dynamical cooling of the systems. Any dynamical cooling due to collisions could pull 
fragments out of resonance and allow them to survive longer to collide with other fragments, 
generating more dust. The timescales and exact evolution of highly eccentric disks, like what 
we would expect from a tidal disruption, is an open question. 

Figure |6] shows the cumulative number distribution of fragments as a function of mass. 
From a close approach distance of 60 to 80 Rwd- Even after a single pass, several frag- 
ments were generated, with an increasing number of smaller fragments as the radius of close 
approach occurred. In all of our simulations, there was a negligible amount of mass lost, 
suggesting that the initial disruption of the asteroid will not show in the WD photosphere 
until dust and material accreted through viscous evolution. 



5. Comparison of the IMMP model to the observed WD population 



We can compare the results of our IMMP model to observed metal polluted and dusty 
white dwarfs. This is possible throu gh the exten s ive co mpil ation of Spitz e r obse rved metal 
enriched white dwarfs presented in iFarihi et al.l (120091 ) and iFarihi et al.l (120101 ) . We took 
these populations and divided them to investigate how each were distributed by WD cooling 
time and the total age of the system. 

The top panel of Figure [7] shows the cooling age for the disk/non-disk populations. 
There appears to be a bimodal population for non-disk systems that cluster around 10 8 yr 
and 10 9 yr, hinting at possibly two different mechanisms, or at least one mechanism that has 
two characteristic timescales after post-main sequence evolution. If we interpret this in light 
of the IMMP model, these two timescales could represent the peak we observed at 30 Myr, 
but scaled to a longer characteristic timescale for perturbation, such as for a more widely 
separated planet. The bimodality could also be a selection effect, since a smaller number of 
hotter WDs have been observed for metal accretion and small number statistics dominate. 
There is not such a clear bimodality for the disk systems. The median cooling age of the 
non-disk systems is 0.9 Gyr while the median cooling age for disk systems is 0.4 Gyr. 

We also calculate the total age of each white dwarf, ttot=^MS+^cooi- The quantity £ms can 
be determined first by inferring an i nitial mass for each WD using an empirical initial-final 
mass function ([Williams et al.l 120091 ) : 
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(M initial - 0.339) 
M final = — (6) 

The main sequence lifetime is then given by t M s = lOMrj^a- Th e bimodality of the 
non-disk systems now disappears, instead the distribution shows a peak around 2.3 Gyr, with 
a high age tail. The disk systems are on average younger, with a peak closer to 1.3 Gyr and 
fewer systems that are older. This may suggest that the underlying mechanism for perturbing 
planetesimals is related to the total age of the planetary system, rather than the time from 
when the mass of the central star changed by a factor of >2, consistent with our IMMP 
model, which depends not only on perturbations, but the total mass of asteroids available 
as a reservoir. The relatively large ages at which disks and pollution appear suggests that 
the mechanism for delivering small planetesimals also must be efficient over Gyr timescales. 
Our simulations show that the IMMP model should be effective at delivering asteroids for 
hundreds of Myr. 

Finally, we can use our IMMP results to infer a distribution of asteroid belt masses 
based on the observed accretion rates for WDs. We assumed that each non-dusty WD 
system represents an "average" accretion rate that can be directly compared to the fit of our 
simulated mass accretion rates. We then determined a scaling factor between our calculated 
accretion rate and that observed to infer a total belt mass using Equation [5j The histogram 
of resulting values is seen in Figure El This may overestimate the masses in a particular 
system if the time between tidal disruptions is longer than the gaseous lifetime of a disk, as 
a disk evolves from higher accretion rates to lower ones. The range in asteroid belt masses 
goes from 4 M$s to 6xl0 5 M$s- The median mass is 820 M$s- This higher mass is may 
not be unusual given that the median progenitor mass of the non-dusty WDs is 2M , and 
the median total age of these systems is 1.5 Gyr younger than our own Solar System. The 
relative youth and greater progenitor mass could correspond to higher belt masses, as could 
differing dynamical histories. This distribution is also biased at higher WD T e ff to the 
largest accretion rates, due to limited sensitivity in the optical to Ca and Mg absorption 
features. More unbiased studies may find a lower median belt mass. We can also compare 
this amount of mass to that i n planetesim als assumed for other models used to explain the 



observed WD accre tion rates. lJural (120081 ) assumed an asteroid belt mass of 10 25 g, and in 
Bonsor et al.l (120111 ). the median mass of their planetesimal belts exterior to a giant planet 
was 10 M e or on the order of 1600 times our assumed belt mass. If our results hold and 
the IMMP model is most efficient at explaining the observations, then it would suggest that 
more massive stars have significantly more massive asteroid belts, or that our asteroid belt 
is depleted relative to the typical planetary system. 
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6. Discussion 

Our simulations suggest a novel and robust way for explaining t he presence of metal 
pollut ed white dwarfs in addition to the mechanism first invoked by iDebes fc Sigurdsson 
(120021 ) and complementary to the scenario suggested by iBonsor et all ( 1201 ll ). While the 
Debes & Sigurdsson and the Bonsor & Wyatt mechanisms require multiple planets and a relic 
icy planetesimal disk to reside within the system, the scenario we present requires only one 
planet roughly the mass of Jupiter as well as a relic asteroid belt similar to the Solar System's. 
Furthermore, since we have looked at only larger asteroids (our sample were incomplete for 
asteroids with R<20 km), there should exist a population of (20/R min ) p (where p is the 
power law value for the asteroidal size distribution and R m i„ is the smallest asteroid size 
that survives post-main sequence evolution) more objects that could be participating in 
tidal disruptions. The recent WISE Mission (IWright et al.ll2010l ) for example will be able 
to better con strain the Solar Sy stem's population of asteroids, as well as find new dusty 
white dwarfs (IDebes et al.ll20f lbl ). We have also shown through tidal disruption simulations 
that significant disruption of smaller asteroids proceeds at radii co mpar a ble to where dusty 
white dwarfs are observed, confirming what had been proposed by I Jural (120031 ) for G 29-38. 
Despite the seeming complexity of such a situation, these scenarios are the best explanations 
for the observations at this time. It is then important to identify which mechanism, if any, 
dominates and whether there are suitable observational tests that could also potentially 
discriminate between planetary instability, the IMMP, or exterior mean motion resonances. 

There are several important implications to our IMMP simulations which will warrant 
further and deeper study. These simulations are useful for understanding the frequency of 
Solar System like asteroid belts as well as their mass, they provide limits to the mass and 
location of planets around dusty/polluted white dwarfs, and they help to explain the relative 
frequency of dusty vs. polluted white dwarfs. 

Firstly, if this mechanism is widespread and more common than perturbations due to ex- 
terior resonances, this directly tests the terrestrial planetesimal population of main sequence 
stars with planetary systems. This is akin to both a determination of %,i an etesimai as well as 
a test of how well populated these regions are. Current observations of main sequence stars 



collis ions or high mass asteroidal belts (e.g. 



Chen et al. 


2006; 


Lisse et al. 


2009; 


Currie et al. 



20 111 , and references therein). Limits on the masses of these types of planetesimal belts will 
be useful for constraining how dusty main sequence stars are in the terrestrial planet form- 
ing regions, a crucial measurem ent to determine h ow successful terrestrial planet imaging 
missions in the future might be (jGuyon et al.l 120061 ). 
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a (All) 

Fig. 1. — The evolution of the libration width about the 2:1 resonance as a function of mass 
lost from the central star (or as a function of planet mass ratio). The libration width grows 
under increasing mass ratio, creating regions that can trap planetesimals that then later are 
perturbed out of the resonance. 
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Fig. 2. — Summary figure of our iV-body simulations of the Solar System asteroid belt under 
the influence of post-main sequence evolution. Asteroids are plotted as a function of e vs. a 
for their initial orbital elements. Overplotted are the libration widths for the 2:1 resonance 
before post-main sequence evolution (solid line) and after the Sun loses half of its mass 
(dashed line). Small asterisks are the asteroids used in our simulations and the squares are 
asteroids that tidally disrupted in at least one of our 17 simulations. The majority of the 
tidal disruptions occur around the edge of the new libration width to the 2:1 resonance after 
the central star's evolution. 
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Fig. 3. — (top) The number of asteroids tidally disrupted as a function of time, normalized 
to the initial number of asteroids in the belt, (bottom) The size distribution of asteroids 
tidally disrupted in our simulations, normalized to the initial number of asteroids in the belt. 
The red dashed line shows the initial distribution of asteroid radii normalized by the total 
number of asteroids used in our simulations, scaled by a factor of 0.0352. 
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Fig. 4. — Observed M vs. T e jj for both dusty (squares) and polluted (diamonds) WDs. 
Asterisks represent the inferred instantaneous accretion rate from all of our simulations 
using Equation HI while the solid curves represent a power-law fit to the simulated accretion 
rates as a function of time with a Solar System asteroid belt (lower curve) and an asteroid 
belt with 820 times the mass of the Solar System asteroid belt (upper curve). 
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Fig. 5. — Snapshots of rubble pile simulations of an asteroid being tidally disrupted by a 
white dwarf. 
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Fig. 6. — Cumulative size distributions of asteroid fragments after initial tidal disruption. 
As the asteroid approaches closer to the WD, more fragments are generated. The distance 
of close approach in WD radii is marked for each simulation. 
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Fig. 7. — (top) Histograms of metal-rich WDs with no disks (solid line) and with disks 
(dashed line) as a function of cooling age (t coo i). (bottom) Histograms of the same types of 
WDs but function of total age, t coo i + t MS . 
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Fig. 8. — Distribution of inferred asteroid belt masses in units of the mass of the Solar 
System's asteroid belt, Mss- Such a distribution could be used to estimate the amount of 
dust present around main sequence stars targeted for terrestrial planet searches. The median 
value of all masses is 820 M S s- 



